#include <iostream>
#include <Eigen/Eigen>
#include <fstream>
using namespace Eigen;

class cal_temperature
{
public:
    cal_temperature(double diff_step, double R_1, double R_2, long int C_in, long int C_wall)
        : dt_(diff_step),
          R_1_(R_1),
          R_2_(R_2),
          C_in_(C_in),
          C_wall_(C_wall),
          P_N_(8000.0)
    {
        A_ << -1 / (R_1_ * C_in_), 1 / (R_1_ * C_in_), 1 / (R_1_ * C_wall_), (-1 / (R_1_ * C_wall_) - 1 / (R_2_ * C_wall_));
    }
    inline void rk45()
    {

        Eigen::Vector2d k_1, k_2, k_3, k_4;
        k_1 = A_ * theta_ + B_ + bias_;
        k_2 = A_ * (theta_ + dt_ / 2 * k_1) + B_ + bias_;
        k_3 = A_ * (theta_ + dt_ / 2 * k_2) + B_ + bias_;
        k_4 = A_ * (theta_ + dt_ * k_3) + B_ + bias_;
        theta_ = theta_ + dt_ * (k_1 + 2 * k_2 + 2 * k_3 + k_4) / 6;
    }
    inline void set_start_theta(double theta_in, double theta_wall)
    {
        theta_ << theta_in, theta_wall;
    }
    inline void reset_B_P()
    {
        B_ << P_N_ / C_in_, 0;
    }
    inline void reset_B_0()
    {
        B_ << 0, 0;
    }
    inline Vector2d get_theta()
    {
        return theta_;
    }
    inline void reset_out_tempurature(double theta_out)
    {
        bias_ << 0, theta_out / (R_2_ * C_wall_);
    }

private:
    double R_1_, R_2_, P_N_;
    long int C_in_, C_wall_;
    double dt_, theta_out;
    Matrix2d A_;
    Vector2d B_, bias_;
    Vector2d theta_;
};

int main()
{
    double R_1 = 0.0012;
    double R_2 = 0.0093;
    long int C_in = 1100000;
    long int C_wall = 186000000;
    double theta_out = 0;
    double dt = 1;
    cal_temperature cal_temp(dt, R_1, R_2, C_in, C_wall);
    cal_temp.reset_B_P();
    cal_temp.reset_out_tempurature(theta_out);
    cal_temp.set_start_theta(0, 10);
    std::ofstream filename;
    filename.open("data.txt");
    for (int i = 0; i < 86400; i++)
    {
        cal_temp.rk45();

        if (cal_temp.get_theta()(0) < 18)
        {
            cal_temp.reset_B_P();
            std::cout << "less than 18";
            std::cout << cal_temp.get_theta() << std::endl;
        }
        else if (cal_temp.get_theta()(0) > 22)
        {
            cal_temp.reset_B_0();
            std::cout << "over 22" << std::endl;
            std::cout << cal_temp.get_theta() << std::endl;
        }
        filename << cal_temp.get_theta()(0) << std::endl;
    }
}